Classical Correlation-Length Exponent in Non-Universal Quantum Phase Transition 

of Diluted Heisenberg Antiferromagnet 



O 

o 
o 

(N 
-i— > 
o 

O 

in 

(N 

a 



C3 



I 

a 

o 

o 



> 

On 

m 
o 

o 
o 

c3 



a 

o 

o 



X 



C. Yasuda 1 , S. Todo 1 , K. Harada 2 , N. Kawashima 3 , S. Miyashita 4 , and H. Takayama 1 
1 Institute for Solid State Physics, University of Tokyo, Kashiwa 277-8581, Japan 
2 Department of Applied Analysis and Complex Dynamical Systems, Kyoto University, Kyoto 606-8501, Japan 
3 Department of Physics, Tokyo Metropolitan University, Tokyo 192-0397, Japan 
4 Department of Applied Physics, University of Tokyo, Tokyo 113-8656, Japan 

(February 1, 2008) 

Critical behavior of the quantum phase transition of a site-diluted Heisenberg antiferromagnet on a 
square lattice is investigated by means of the quantum Monte Carlo simulation with the continuous- 
imaginary-time loop algorithm. Although the staggered spin correlation function decays in a power 
law with the exponent definitely depending on the spin size S, the correlation-length exponent is 
classical, i.e., v = 4/3. This implies that the length scale characterizing the non-universal quantum 
phase transition is nothing but the mean size of connected spin clusters. 
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Ground-state phase transitions in two-dimensional 
(2D) diluted quantum Heisenberg antifcrromagnets 
(HAF's) have attracted much interests because they are 
caused by the coexistence of quantum fluctuations and 
randomness • Numerical works |^|,^| and theoretical 
works fl(],[ll| have given various estimates of the critical 
concentration of the system with S — 1/2. All of them 
are above the purely-geometrical percolation threshold 
on a square lattice, p c i=0- 5927460(5) [Q. This suggests 
that the phase transition could be seriously affected by 
quantum fluctuations. 

Recently, Kato et al. |l3| have investigated the di- 
luted 2D HAF with 5 = 1/2, 1, 3/2, and 2 by means 
of the quantum Monte Carlo (QMC) method with the 
continuous- imaginary-time loop algorithm Jl4|-4j7| . Their 
conclusion is qualitatively different from the above men- 
tioned results: the critical concentration coincides with 
Pd and does not depend on 5. The coincidence of the 
critical concentration with p c \ has also been reported on 
the bond-diluted HAF @. 

The critical exponents of the phase transition at p c \ 
have also been estimated. Kato et al. have obtained 
the critical exponent (3 of the zero-temperature staggered 
magnetization. Interestingly, the value of (3 is different 
from that of the classical (5 = oo) exponent and depends 
on 5. They have also estimated other critical exponents 
by the finite-size scaling (FSS) analysis assuming the fol- 
lowing form for the static staggered structure factor at 
zero temperature: 

5 s (L,0,p)~L*5 8 (L 1 /"(p-p cl )) , (1) 

where L is the system size and 

S s (L,T,p) = J-^e^W-^SfS;) (2) 

with k = (it, 7t) and d = 2. The bracket (• ■ •) in Eq. 
(0) denotes both the thermal and random averages. The 



scaling function S s (x) in Eq. ([!]) has the following asymp- 
totic form: 



S s (x) 



x 2 P for i>1 
Id""* for x < -1 , 



(3) 



where the exponent (3 is related to 'F and v by the scaling 
relation 



2/3 = (d- V)v 



(4) 



They have estimated the critical exponent \P by the FSS 
analysis exactly at p = p c \ and have found that it depends 
onS. 

Kato et al. have attributed the 5-dependence in ex- 
ponent (3 or 'F to quantum fluctuations, while the length 
exponent v is assumed to be given by the classical one, 
v c \ = 4/3 [fl9|], which governs a power-law divergence 
of the mean size of connected spin clusters as X(p) oc 
\p — Pci\~ VcX - Namely, they have assumed that the stag- 
gered spin correlation function between two spins in a 
cluster is described by the scaling expression 



C(i,j;p)~r7fC(r i , j /\(p)) 



(5) 



with G(x) ~ const, at i < 1. The power-law decay of 
the correlation function is due to quantum fluctuations 
and its S'-dependent exponent a is related to 'F by 



* = 2D - d - 



(6) 



where D is the fractal dimension (91/48 for d = 2). Kato 
et al. have checked this scenario simply by evaluating v 
through Eq. (Q) using (3 and \P obtained by their simu- 
lation. 

In the present paper, we perform the FSS analysis on 
S s (L,0,p) of systems with 5=1/2 and 1 more system- 
atically by carrying out the QMC simulation at various 
concentrations not only at p c \. It is found, as we see be- 
low, that our QMC data are well described by the FSS 
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form of Eq. (@) in a whole range of p studied, which 
includes the asymptotic ranges x 3> 1 and x <C —1 in 
Eq. ([}|). Within numerical accuracy of our analysis, ex- 
ponent 'J turns out to definitely depend on the spin size 
S, while exponent v does not. The latter S'-independent 
value of v coincides with the classical one (=4/3) within 
the error bar. These results further support the above 
mentioned scenario of the quantum phase transition in 
the diluted 2D HAF, particularly, the ansatz that the 
length scale characterizing the transition is nothing but 
the mean size of connected spin clusters. 

The system we study is the site-diluted HAF on a 
square lattice described by the Hamiltonian 



(hj) 



(7) 



where J(> 0) is the antiferromagnetic coupling constant, 
53 (ij) denotes the summation over all nearest-neighbor 
pairs and Si is the quantum spin operator at site i. 
The quenched magnetic occupation factors {e^} indepen- 
dently take 1 or with probability p and 1 — p, respec- 
tively. We simulate L x L square lattices with the pe- 
riodic boundary condition by means of the same QMC 
method with the continuous-imaginary-time loop algo- 
rithm as that adopted by Kato et al. An improved esti- 
mator is used to calculate the static staggered structure 
factor. At each parameter set (L,T,p) physical quanti- 
ties of interest are averaged over 1000 - 3000 samples. 
At each sample, 10 3 - 10 4 Monte Carlo steps (MCS) are 
spent for measurement after 500 - 10 3 MCS for thermal- 
ization. 

In a finite system, S a (L,T,p) converges rapidly to its 
zero-temperature values at temperatures lower than the 
gap due to the finiteness of the system. The satura- 
tion temperature turns out to be smaller for smaller 
\p — p c i\- In the present work S s (L,0,p) close to p c \ is 
approximated by S s (L,T,p) at the following tempera- 
tures: for the S = 1/2 system T = 0.002 J for L = 24 
and T = 0.001J for L = 32, 40, and 48, and for the S = 1 
system T = 0.01J for L = 24 and T = 0.005 J for L = 32, 
40, and 48. 

Let us first examine the FSS analysis of S s (L,0,p) at 
p close to p c \. The results of the S — 1/2 (S — 1) system 
at 0.580 < p < 0.605 (0.585 < p < 0.600) are shown in 
Fig. 1(a) (1(b)). The error bars in the figures represent 
the standard deviation. In the FSS fit the critical con- 
centration is set to be p c \ (~ 0.5927460), and the scaling 
function S s (x) is approximated by a polynomial of order 
2. As seen in Fig. 1(a), the QMC data for S = 1/2 are 
well scaled with = 1.192 and v — 1.23. The statistical 
accuracy of the fit is shown in Fig. 2, where we draw 
the confidence region within which the true values of W 
and v fall with probability 68.3% (1-er), 95.3% (2-cr), or 
99.7% (3-<t). Similarly, the data for S — 1 are well scaled 
with \I> = 1.555 and v = 1.09 as seen in Figs. 1(b) and 
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FIG. 1. Scaling plot of S s (L,0,p) for (a) S = 1/2 and (b) 
S = 1. The dashed line represents S s (x) which is approxi- 
mated by a polynomial of order 2. 



2. From these results we can conclude that exponent VP 
definitely depends on the spin size S. The exponents v 
of S = 1/2 and 1, on the other hand, coincide with each 
other and with its classical value (= 4/3) within the nu- 
merical accuracy of the present analysis. The obtained 
values of v, a, and (3 are summarized in Table I. The 
exponents a and (3 are calculated by the scaling relations, 
Eqs. (|) and (§). 

The value of \& for S — 1/2 obtained above is consis- 
tent with the one (= 1.17(6)) that has been estimated 
by Kato et al. from S s (L,0,p) at p = p c \. In their anal- 
ysis S B (L, 0,p c \) is approximated by S s (L,T,p c {) at low 
temperatures where its T-dependence becomes not dis- 
cernible within the error bars. Kato et al. have also 
performed the FSS analysis making use of the data at all 
temperatures they have simulated. This analysis yields 
^ = 1.27(2) which differs distinctly from the present re- 
sult listed in Table I. This discrepancy may be due to 
the systematic error in our estimates where S s (L,0,p) is 
approximated by S S (L, T,p) at a small but finite temper- 
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FIG. 2. Confidence region of $ and v. The percentages 
68.3% (1-cr), 95.4% (2-cr), and 99.7% (3-cr) represent the prob- 
ability that the true parameter values fall within the con- 
fidence regions. The classical value of ^ is 1.79167. The 
vertical dashed line indicates the classical value of v, 4/3. 

ature as described before. It should be remedied when 
the finite-temperature FSS analysis is carried out as done 
by Kato et al. 

From Fig. 2 one sees that the accuracy of v is sig- 
nificantly poorer than that of which is interpreted as 
follows. The value of iff can be essentially extracted solely 
from S S (L, 0,p c i), while to evaluate v we need S S (L, 0,p) 
at p other than p c \, or SS s (L,0,p)/Sp at p = p c \. Nat- 
urally, the statistical error of v is expected to be larger 
than that of \& . Concerned with the systematic error 
in v arising from finite-temperature corrections, on the 
other hand, the saturation temperature of S S (L, T,p) be- 
comes smaller for smaller \p— p c \\, i.e., the effect of finite- 
temperature corrections to S s (L,T,p) becomes maximal 
exactly at p = p c \. This implies that the leading order 
finite-temperature correction to 5S S {L, 0,p)/Sp at p = p c \ 
disappears, and so the systematic error in v is much 
smaller than that in iff. 

Next let us discuss S S (L, 0,p) in the full range of p we 
have simulated, i.e., 0.2 < p < 1, which is much wider 
than that in Fig. 1. As shown in Figs. 3(a) and (b) re- 
spectively for 5 = 1/2 and 1, all QMC data turn out to 
lie on a universal curve when L~^S S (L, 0,p) are plotted 
against L 1 ^^— p c i\ by using the exponents iff and v listed 

TABLE I. Summary of critical exponents iff, u, a, and f3 
for S — 1/2, 1 and the classical case (S = oo). The values 
of ^ and v are obtained by the FSS analysis shown in Figs. 
1(a) and (b), and a is calculated as a = 2D — d — if 1 , and P 
as /3 = (d-ty)v/2. 
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FIG. 3. The double-logarithmic plot of L' 9 S S (L, 0,p) 
against £ V "|p - p c i| for (a) 5 = 1/2 and (b) 5=1. Dashed 
lines represent ax 213 and 6|a;| - "*. 



in Table I. For L 1 ^^ — p c \\ > 1, the data points merge 
to the dashed line, which represents ax 2 ^ for p > p c \ and 
&|a:| -I/ * for p < p c \, where /3, v, and iff are those listed in 
Table I, and a and b are arbitrary constants adjusted to 
fix the position of the asymptotic lines. At much larger 
l// v \p—p c \\ the scaling fit becomes deteriorated, indicat- 
ing that the corresponding p is out of the scaling region. 
The results shown in Fig. 3 also support the scenario 
due to Kato et al. on the non-universal quantum phase 
transition. 

In summary, in order to establish nature of the non- 
universal quantum phase transition of 2D site-diluted 
HAF's for S = 1/2 and 1, we have performed the QMC 
simulation in a relatively wider region of concentration 
than that of Kato et al. and have estimated the critical 
exponent v more systematically. We have observed that 
the static staggered structure factor is well described by 
the scaling form of Eq. (|l|) with Eq. (|J). In particular, 
the exponent v is confirmed to coincide with the classical 
one. These results support the following arguments by 
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Kato et al: 1) there exists no other macroscopic char- 
acteristic lengths than X(p), the mean size of connected 
spin clusters at concentration p and 2) the staggered spin 
correlation function between two sites on a fractal cluster 
decays in a power law as C(i,j;p) ~ r 7?> where a de- 
pends on the strength of quantum fluctuations specified 
by the spin size S. 

Most of numerical calculations for the present work 
have been performed on the CP-PACS at University 
of Tsukuba, Hitachi SR-2201 at Supercomputer Center, 
University of Tokyo, and SGI2800 at Institute for Solid 
State Physics, University of Tokyo. The present work is 
supported by the "Large-scale Numerical Simulation Pro- 
gram" of Center for Computational Physics, University 
of Tsukuba, and also by the "Research for the Future 
Program" (JSPS-RFTF97P01103) of Japan Society for 
the Promotion of Science. 
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